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JU INTRODUCTION 


Recent research at the NPS in atmospheric global predic- 
tion has been primarily concerned with the five-level global 
primitive equation model initially designed and programmed 
by Dr. F.J. Winninghoff and further developed by Elias (1973), 
Mihok (1974), McCollough (1974) and Maher (1974). Their 
numerical experiments proved valuable but pointed out the 
inflexibility of the model for research purposes. 

Arakawa and Mintz (1974) described an atmospheric global 
prediction model with a new finite difference formulation 
and greater vertical resolution. The model employs a new 
horizontal distribution of variables to improve the geo- 
strophic adjustment process. The purpose of this research 
was to program and check out a modified version of the new 
Arakawa scheme. The fundamental approach was to develop a 
flexible program for use in a number of projects including 
modeling of tropical circulation. Eventually the model's 
performance can be compared with the performance of FNWC's 
global model. 

The model was initially developed with two levels and a 
variable grid to aid in debugging. The atmosphere was 
treated as adiabatic and frictionless and the water vapor 
continuity equation was not considered. The initial condi- 
tions for all experiments were analytically derived in the 


same manner as Maher (1974). The advantages of analytic 
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initial conditions were significant in reducing computer time 


during check out of the model. 
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II. MODEL DESCRIPTION 


The differential equations are essentially the same as 
those described by Arakawa (1972). The integration scheme 
was carried out on a staggered, spherical, sigma coordinate 
system. For an adiabatic, frictionless model the primitive 
equations form a closed set. No sink or source terms were 


considered in the course of development or check out. 


A. VERTICAL COORDINATE SYSTEM 

The model uses the non-dimensional sigma coordinate sys- 
tem as described in Haltiner (1971). The two levels divide 
the troposphere in half and the tropopause is assumed iso- 


baric. The sigma coordinate is defined as 


b š 
fo 7 (2.1) 


where p is the pressure, p+ the constant tropopause pressure 
and T is the terrain pressure. The terrain pressure is fur- 


ther defined as 


qom Ps a P+ 1 (2.2) 
where D. is the surface pressure. It follows from equation 
(2:1) that 

G = O at p = Py > 

(2.3) 

G = 1 at p = Ps ° 


which are the vertical boundaries of the coordinate system. 
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n 





Figure 1. The sigma (o) coordinate system as used in the 


mode1. Py is 200 mb. 
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The boundary condition at 0=1 is 0=0 and at 0=0 it is as- 


sumed that 0=0 where o is defined as do at. 


B. PRIMITIVE EQUATIONS 

The continuous form of the primitive equations is de- 
scribed in the orthogonal curvilinear coordinates & and n 
where & = à (longitude) and n = Q (latitude). An area ele- 


ment is Z Ag An and the lengths of the sides are = and an ; 


where 


I 
= = a. 
n 


E 


= a cos h and 


The difference element is illustrated below. 


17 + A?) | | 


; 


Š AC £ + AŠ 
m 
The mass continuity equation is 
3 , T 3 , u 2 3 TO. 
3t "mn? + 3g U 8) * anm) * Soma) ° 9: ED 





The &-component of the horizontal equation of motion is 





d T O -TU O `, TV 3 ,TG 
a E praece p mmu) K 9c n 
f o 1 o 1 T- oğ Ə Tr T 
E — + — — -— — — — | —— — = —əƏ 
Z n S R E na p (2.9) 


The n-component of the horizontal equation of motion is 


d T 9 ,TU 9 ,TV 3 ,TO 
pica | «(mm >) © damn”? 
f o 1 OBEN T, 96 Om M T 
m ` v7: n “n OMS T N S ^ mn a Zn 


The first law of thermodynamics is written in the following 


manner: 

C, moon * Q, CD 
where 

w = E = 10 + Sz - V*V)m (2.8) 


and the flux form is 


9 > Den: _ f 
ED + C (TV > + TUI = T(wa + Q). (2.9) 
Using the relationship for potential temperature, 
K 
T = 9) (2.10) 
0 


the final form in curvilinear coordinates is 
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Pour 3 (mu 3 (mv pa o nS 
m Cp? T SES E B 3n m Fe = Gm on > 


A DN. u T _ von ES 
= Toal J+ BE t m 5 T mx Q . | (2.11) 


The equation of state is 


RT 
= — 2419 
ol D ( ) 


and the hydrostatic relationship is 
Gd = -maóo . 22193) 


A complete set of symbols for the above equation may be 


found in the front of this report. 


C. VERTICAL INDEX 

The index k is used to identify the vertical levels. At 
the upper boundary k=0, P=P, and at the lower boundary 
k=K+1 and P=P.- The variables V and T are carried at the 
odd levels while no is carried at the even levels. The fol- 
lowing definitions are required to complete the vertical 


system: 


Ao 


E (2.14) 


k ` S k+1 k-1 ' 


K 
> ' 
yii ^9, | (2.15) 


where the summation is over odd k (see figure 2). 
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INDEX COMPUTED SIGMA 











k VARIABLES | 
0 TOO oczo 
1 YT 
2 TO O = 0, 
3 V T 
k-2 VT 
k-1 TO Co 
k L. ç 
k41 mo g = + 





k+2 i s WT oO 











K-2 VT 
K-1 TO o=0,_, 
K VT 
K +1 TO=0 Om 1 
VERTICAL INDEX 
Figure 2. In the above figure k is a variable vertical 


index, sigma (o) is the dimensionless vertical coordinate, 

- * . ° ° . 

V is the horizontal vector velocity, T is the terrain pres- 
sure, G is the vertical velocity, T is the temperature and 


® is the geopotential. 
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D. HORIZONTAL GRID AND DISTRIBUTION OF VARIABLES 

The seleetion of the horizontal distribution of variables 
is: a function of two distinct processes. The first is proper 
simulation of the geostrophic adjustment process and the 
second is proper simulation of slowly changing quasi-geo- 
strophic motion after it has been established by geostrophic 
adjustment. As shown by Winninghoff (1968) geostrophic ad- 
justment depends on how the variables are distributed over 
the grid points. Winninghoff used the following equations, 
which are the simplest ones in which geostrophic adjustment 
can take place, to demonstrate five possibilities for the 


placement of the dependent variables: 


du - oh 


dt fy G O, (2.16) 
dv = oh 

== f — 

dt fu + Say O, (2 17) 
dh dU OV. _ 

a Az * dy” S G 


The above equations represent an incompressible, homogenous, 
non-viscous, hydrostatic, rotational fluid with a flat bot- 
tom and a free surface, where u and v are the velocity com- 
ponents, h is the depth of the fluid, f is a constant 
Coriolis parameter, t is the time, x and y are the horizon- 
tal coordinates and g is gravity. The five possible distri- 
butions of the dependent variables h, v and u are on in 
figure 3. Scheme B was used in the Mintz-Arakawa two-level 


global circulation model (Langiois and Kwok, 1969) and 
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i - 1/2 í i * 1/2 


Figure 3. The above figure represents the placement of the 
dependent variables h, v and u where h is the depth of the 
fluid, u and v are the velocity components and d is the grid 


distance (for the simplest case of geostrophic adjustment). 
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scheme E is used in the global model now under further de- 
velopment at FNWC (Maher, 1974). 

In the one-dimensional case it was shown that both 
scheme B and C adequately simulated geostrophic adjustment. 
In the two-dimensional case, however, scheme C was shown 
most satisfactory to simulate geostrophic adjustment, except 
where Ap/d is less than or close to one. The quantity An is 
the Rossby radius of deformation and d is the grid distance 


(Arakawa and Mintz, 1974). The Rossby radius of deformation 


is 


where H is the mean value of h, g is the acceleration of 
gravity and f is the constant Coriolis parameter. The con- 
dition Ap /d Eu cancabuormalocase, therefore, the horizon- 


tal distribution of variables is based on scheme C. (See 


figure 4.) 


E. TIME DIFFERENCING 
The time differencing is carried out in thirty minute 
Sequences. The initial step in each sequence is a two part 


Matsuno scheme represented by the following notation: 


t 
K n 

EE At = (Forward), (2.19) 
t+1 t 338) | 

F = F° + At 3t (Backward) . (2.20) 
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"—] | x | | 


i-1 i-1/2 i i+1/2 id 


Figure 4. The above figure represents the horizontal distri- 
bution of dependent variables with polar modification. Pin 
DO system erepresents the north pole, 1 represents 
the south pole and i represents a meridian. T represents 

the variables T, ® and wm carried at "w-points." u and v are 


the horizontal components of velocity. 
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F is a vector representing the dependent variables. The 
superscript t represents the time step, the superscript * 
represents the results of an intermediate step and At is the 
time interval of a single step. The Eins time steps in 


a sequence use the leapfrog scheme as follows: 


t 
Ug UL Ont — (2.21) 


At the end of each sequence provision is made for the calcu- 
lation of the source terms. The present model doesn't in- 


corporate source term calculations. 


F.  AVERGING THE PRESSURE GRADIENT AND ZONAL MASS FLUX NEAR 

THE POLES 

A problem that may arise at higher latitudes is computa- 
tional instability. Computational instability is a result 
of the convergence of the meridians to the poles causing 
greatly reduced grid distances along a latitude circle. 
Therefore, some technique must be employed to eliminate the 
stability problem. One technique involves reducing At in 
higher latitudes but it would impose serious programming 
difficulties and require more computer time. Another tech- 
nique would be to reduce the number of grid points in a 
latitude circle at higher latitudes, however, it also would 
impose programming difficulties. A third technique is zonal 
smoothing which was developed by Arakawa and Mintz (1974). 
Smoothing preserves the integrity of the horizontal grid 


while maintaining a constant time step. Zonal smoothing is 
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used in this edge und also by FNWC to eiiminate computa- 
tional instability at high latitudes. 

Zonal smoothing involves expanding the pressure gradient 
and zonal mass flux into a Fourier series and reducing the 


amplitude of each wave component by a factor 


= D cos o 
9 7 Gt sin (md) CIL 
where 
S = Stability Coefficient, 
D = Grid distance at the equator, 
Co = Phase speed of the fastest gravity wave, 


At = Time step, 
m, = Wave number, 


d = Grid distance in degrees. 


If the value of S is greater than one, smoothing is not re- 
quired (Arakawa and Mintz, 1974). As pointed out by Arakawa 
the smoothing operation does. not smooth the fields of vari- 
ables, because it is simply a generation of multiple point 
difference quotients. The bar operators in Chapter III in- 


dicate zonal smoothing. 


G. PROGRAM FORMAT 

The program format was patterned after the Mintz-Arakawa 
two-level global model as described by Gates, et al. (1971) 
and Langlois and Kwok (1969) and modified to include scheme 
C distribution of variables and a Matsuno-leapfrog time in- 


tegration scheme. Figure 5 gives a simplified flow diagram 
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for a five-day forecast. The main program controls the 
overall time period of the forecast as well as input and 
output. Time integration for each 30-minute sequence is 
controlled by subroutine Timestep. The forcing functions 

are calculated in Comp 1 and Comp 2 with zonal smoothing ap- 
plied as necessary. The remaining subroutines are peripheral 
to the main flow of the model. 

A main consideration in the program was flexibility, 
therefore, the vertical structure along with the horizontal 
grid are variable. The model can be used to simulate sever- 
al different horizontal and vertical structures including 


the FNWC's 5-level global model. 
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III. FINITE DIFFERENCE EQUATION 


The finite difference equations given in this chapter 
were developed by Arakawa and Mintz (1974). The advection 
terms of the equation of motion were modified to eliminate 
the diagonal flux calculations. The short-term performance 
should not be affected by elimination of the diagonal flux 
terms and considerable computer time was saved. 

The notation used in this chapter represents the specif- 
ic variable centered in the ij index system. Figure 6 con- 
tains an example of 'm-centered'" and "u-centered' notation. 
T-centered notation is used for the continuity equation and 
the thermodynamic equation in which the dependent variable 
is carried at a m-point. U-centered notation is used for 
the £-component of the equation of motion and v-centered 
notation is used for the n-component of the equation of 
motion. Refer to the front of the report for a complete 


list of symbols. 


A. CONTINUITY EQUATION 
The following form is used for the continuity equation 


given in equation (2.4): 


oli. 
J ak Em uk ERE 

9 t edo] "Ket, J Si, 343 St jud 

Hc we Ski) - 9. (348) 
Mak 1:3 1,J 
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^ T-CENTERED 
NOTATION 





i u-CENTERED 
| | NOTATION 


Figure 6. An example of the notation used to describe the finite dif- 


ference equations. The continuity equation is described on a "m-centered" 


grid in which the F and G symbols are flux calculations in their re- 


spective directions. The "u-centered" grid is an example used to de- 


: : U U 
scribe the &-component of the equation of motion where F and g are 


also flux calculations. 
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(bar indicates zonal smoothing) 








K 
: = E Er m A : 
Oi, itd T 543 Ex 76 2, 
= ALAN 
Il : T (===> 3 
1,J 1,J mn s 
S =. .c 
drd O E 
The tendency equation is 
ONES o Ie 
o eS M We ee, 
ot rou. J ls QN 1, jee x 
The vertical motion equation is 
Sel K, ji E. m E ye 
359) 5 1 142,0 ed 1. ji 103 2 
Bui. 
_ štit dh (5.5) 
9t 
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B. PRESSURE GRADIENT Term AND HYDROSTATIC EQUATION 


s Hydrostatic Equation 


The geopotential for a level K is given by 


K RT. . 
K S t k i K 
ES ?,)99 x51 "GG? Tk g 
P1, j 
K-2, R K2 K k K 
MEM uc EEN B | (3.6) 
k=1 pec j p : É 
| IE M 13 


where 9^ is the geopotential at the surface, 





^k*1 . in 8° - In 0 
E 
5Kt2 ~ GK 


Given n the gecpotential for the remaining levels is found 


using the following relationship: 


k+2 K k K 
k k+2 Ak+1 
DM o Ó CE >) A 91 


m (3.7) 
asd 1,5 P Po “ij e 





2.  n-component of the Pressure Gradient Force 
The pressure gradient force in the meridional equa- 


tion of motion as given in equation (2.6) is 


E AL 4 k _ k 
nr BLA jag * 71,3742 053,344 T 91, 3-2? 
+ ((noa)* + (noo) . ¿Mr - om )] (3.8) 
Oi tk i,j-2/  i,j*à i,j-à^! ' j 
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where 


| ors = T G 
i, jtó aus p 


3. E-component of the Pressure Gradient Force 


The pressure gradient force in the zonal equation of 


motion as given in equation (2.5) is 


AT a k EN K 
EE Ds 6s 0o3 7 757 
k k 
+ ((700):,4 E (moa), mr E Ti 2M : (3.9) 


Oe FO ENTUM ELUXES ` 
The zonal momentum flux as given in equation (2.5) is 
9 U k 1 U k 
3t 04,4 1,32 * 8L £444, 30444,3 * 1,3? 
u k u k 
_ š ET o : * E. s B T ES - 
E PEE, ES Yi, 3? 
u k ioe A K+2 k 
_ ° M = se E + — > n uv or rem 
Bi ja, al? ] Ach LE Cu; j u; j) 
A Y, (3.10) 
wJ 2x 2-5) 
where 
u RE 
MS y alita, j M H; y j) , 
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pe = NA + op* m p* I ! 
DA) “` j+£, j+1l i+3,j Ban , (3.11) 
u cod * * k * T * 
Bi jag ^ (G uy Mt Ct GC, uu (83.12) 
where 
p“ - d E + F oco 
A A? c (3.13) 
G. . = 4G + G 3.14 
i,j = 206; jes 7 Sa, 5-2? (3.14) 


The meridional momentum flux is Similiar to equation (3.10) 


with u replaced by v. (les ay, ee 
form: 

V En = T 
Fira, 3? ia, 39 7 71,548 

V TP E * 
Ej j+ ^ Gir jr jr 

V LU E 

a AR 

e Y ° ° 

S. . = %(S. . + S. . 

i,j ^ $01,344 i,j-2)2 

* 
F is defined by 
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and SV take the following 


* 


Ee eus c c a 

G? ) (3.16) 
i-1,j*2$ ; : 

(3.17) 

(3.18) 


* 
equation (3.13) and G by equation (3.14). 
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E RIRI LG FORCE 


A variable e F is defined at m-points by 


5 


= A&An e NES 
NEED m ie 


E ] 
Jz 


JTA 2 


where L. 7 is the Coriolis parameter. 


, 


The zonal-component of the Coriolis force as given in 


equation (2.5) is 


k k 
4 + 
LEE l 1253 EMI Vi+z, 3-3) 
k k 
E tov. I 
IE a io] ° d) 


k k 
-4In. . GR: u + u 
al Ju jor llo.) i-k,j-2) 
k E k 
a . + . l ] l 
T Ci j+2 (Uitg, j+ E ] (3 20) 


E. THERMODYNAMIC ENERGY EQUATION 
The thermodynamic equation corresponding to equation 


(2.11) may be written as 
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POLAR MODIFICATION 


The poles of a spherical coordinate system are singular 


points and polar velocity components cannot be defined. 


Terrain pressure (T) at the poles can change as a result of 


meridional flux at all points P-Z and P+% where v is carried. 
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ITI HT equation is modified at the poles by omitting 
all undefined flux terms. Each pole is treated as a series 
of points in order to facilitate programming (see figure 7). 
Each index i,P represents the shaded area B in figure 8. 
The continuity equation is integrated at each point repre- 
Senting the pole and then averaged to determine a polar val- 
ue. The thermodynamic equation and the vertical velocity 
are treated in a similar manner. The advective terms in the 
equation of motion are given special treatment. The remain- 
ing terms in the equation of motion are unchanged except in 
the case where they are undefined and omitted. 


1. The Advective Terms in the €-component of the Equa- 


tion of Motion 

The polar modification of the ¿-component of the 
equation of motion is shown in figure 9. It was modified 
from the Arakawa scheme by eliminating the diagonal flux 
terms. The following form is used: 


*u k 
s 
| Ma C p.34 * 9444 p) 
i,P-1 : ; 


*u 


k u 
pa SA, pay * a.p) -€ 


k 
o Ap E pee 


1 U,k+1 


EON k+2 k G: k-2 
TE eps ee uu p» Sus scm EA f: 
(3.22) 
where 
*u EE * T * 
E 90S p, * Pig p»: 
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Figure 8. Bach index 1,P in the continuity equation is 


represented by the shaded area. m at the poles 
can only change as a result of G. The thermo- 


dynamic equation and vertical velocity are treated 
in the same manner. 
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Figure 9. 





The polar modification of the u equation of motion. 
* 

F Y and g? are flux terms. The shaded portion 

represents the area associated with each variable 


u, 
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gu given by equation (3.12) 


u = 
ee sp) > 
S = 4(Š + § ) 
oes eee eed 


2. The Advective Terms in the n-component of the Equa- 


tion of Motion 

The polar modification of the n-component is repre- 
sented in figure 10. It was modified from the Arakawa 
scheme in a Similar manner as equation (3.22). The follow- 


ing form is used: 


9 (II Vv ) 1 * yv r: k 
* v / k V k 
Noe Eos * Vi oq p.;/ By OUR ose po) 
INL v ktl, kt2 k š «NECK k-2 
—— es Rn 3 + = > n 
+ n z[ S NEP, PRESS Oa eo vi pao 
(3.23) 
where 
Dy E at * 
z t 
n p 5 pn a 


gY given by equation (3.16), 


= 1 + 
Hi: poz I, pt Hy py) » 
e V - x ⁄ o 
TE 
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Figure 10. The polar modification of the v equation of 
* * 
motion. F Y and g V are flux terms. The shaded 
portion represents the area associated with each 


variable v. 





IV. ANALYTIC INITIALIZATION 


The initial conditions were similar to those used by 
Maher (1974). Analytic initial conditions offer several ad- 
vantages in the initial evaluation of a model; they simplify 
the otherwise difficult task of balancing and interpolating 
initial conditions from constant pressure surfaces to sigma 
(o) surfaces, they allow the angular phase speed to be esti- 
mated from a non-divergent model and they allow the simula- 
tion of atmospheric states that otherwise would be difficult 


to simulate. 


A.  ANALYTIC BALANCING 

The initial terrain pressure and velocity ielas were 
obtained from the stream function solution to the linearized, 
non-divergent vorticity equation (Haurwitz, 1940). The com- 
plete vorticity equation was later solved by Neamtan (1946). 
The stream function, y, may be written as follows: 

lr 
y 7 A sin(m,A - vt) sin à cos ~ $ - Ba? sin dh. (4.1) 


where A and B are constants, v is the angular velocity, m, 


is the wave number and a is the radius of the earth. The 


angular phase speed is given by 


EXEC NOI) = 2 20 


n- ^ P^ NOD C NOM) ' po 


where v/m, is the angular phase speed, N-m, *1 and Q is the 


angular velocity of the earth.  Harmonic waves defined by 
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the stream function will move with a constant angular veloc- 
ity without changing shape assuming a baratropic, non-diver- 
gent atmosphere (Haurwitz, 1940). However, the equations 
used in this model are for a divergent atmosphere and will 
give a smaller contribution from the second term in the equa- 
tion (4.2), especially with small wave numbers. 

Equation (4.1) is used as the forcing function to obtain 
the geopotential from the non-linear balance equation 


t 


(Phillips, 1959). The geopotential perturbation, % , may be 


written 
ie a*A(>) + a*B(d) sin mA + a*C(o)(2 sin*m,A re 
where 
2m 
B à A ^ k , 
A($) = 5(22 + B) cost + az) cos `$ 
2m, 
2 2E _ _ 
[ (m, *1) COS Q + (2m) a 2) LT : (4.4) 
A 
2(8+B )=7 m. à 
B(¢ ) Ma) °° Q [(m, + 2m, * 2) 
k 
=: (m, *1)* cost bl. (4.5) 
A 2 2m. 5 
C(9) = $(z27) cos $ [(m,+1) cos’$ - (m, +2)] . (4.6) 
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B. ANALYTIC WINDS 


The wind components were obtained as follows: 


Pee ot) | 
lE er E (4.7) 


ec (4.8) 


a cos $ ddA 


Performing the operations indicated above the u and v com- 


ponents may be written 


m. +1 


ag t. _ k 
u = a LA sin(m,A vt) Sos ó - m, A 
m, -1 
sin(m À - vt) cos @ sin-¢ — Ba- cos ¢]_, (4.9) 
Y a T ne sus 4 Rowe ó cos (m À — vi)] (2110) 
SU f k f ; 


C. TEMPERATURE 
The temperature field was derived in accordance with the 
NACA standard atmosphere as follows (Haltiner and Martin, 


1957): 


T(°K) = 288 - 0.00652, (4.11) 


02217970723 


Z (meters) = 44308 [1 - (43-55) 


I 2) 


where CE TOGO meters, 
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V. MODEL PERFORMANCE 


The evaluation of the model was a two-step process. The 
first step involved the use of a coarse grid and was essen- 
tially a scheme to minimize computer time while debugging. 
The second step involved the use of a finer grid for better 
resolution. The majority of the experiments presented in 
this chapter were run with the fine grid. 

The coarse grid consisted of 16 points N-S which gave a 
12 degree separation between grid points. Since strictly 
analytic initial conditions were used it was only necessary 
to integrate over one wave using cyclic continuity and there- 
fore a 10-point E-W grid was adequate. The E-W grid dis- 
tance was a function of wave number. Wave numbers 4, 8 and 
12 were used in the experiments. Wave number 4 gave an E-W 
grid spacing of 9 degrees, wave number 8 gave a grid spacing 
of 4.5 degrees and wave number 12 gave a grid spacing of 3 
degrees. The fine grid used 46-points N-S which gave a grid 
spacing of 4 degrees. The E-W grid was similar to the 
coarse grid. 

In all cases the model was run with two levels, a flat 
earth and no source or sink terms. The time step (At) was 
six minutes. Fourier analysis of the surface pressure field 
was used to compute phase speed and wave amplitude in all 
but experiment III in which a Fourier analysis of the wind 


field was computed. Calcomp charts were plotted at 12 hour 
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intervals and consisted of one wave of the surface pressure 
field over the complete N-S grid. The phase speed and wave 
amplitude were initialized through the constants A and B 
given in equation (4.1), where A is the amplitude of the dis- 
turbance and B is the amplitude of the mean flow. 

Experiment I. Wave numbers 4, 8 and 12 were used to de- 
termine the phase speed over the coarse, 16-point grid. The 
phase speed was set at 10 degrees per day by adjusting the 
amplitude of the mean flow (B). The graphs of phase angle 
(degrees longitude) vs latitude are shown in figures 11, 12 
and 13. The observed phase speed approached 20 degrees over 
48 hours for all three cases. 

Fopoepiment IIT. This and all subsequent experiments were 
run over the 46-point N-S grid, External gravity waves were 
simulated with wave number 4 by setting the mean flow, the 
u and v components of the wind and the Coriolis terms to 
zero. The theoretical period of the gravity wave was calcu- 


lated using the following approximation: 


RASO OO ) 
E = ons.” C5.» 
p 
where 
K = YN(N+1) 
p a 


C is the zonal phase speed (set at 300 m/sec), E is the 
two-dimensional wave number, N=m, +1, a is the radius of the 


earth and E is the period. The above approximation gives a 
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period on the order of 7 hours for wave number 4. The ob- 
served oscillations of the surface pressure field as a func- 
tion of time are given in figure 14. The period is on the 
order of 6-7 hours. 

Experiment III. In this experiment pure advection was 
Simulated with wave number 4 by removing the Coriolis terms, 
the pressure gradient terms and the vertical advection 
terms. The surface pressure and the temperature were kept 
constant. The u and v components of the wind were Fourier 
analyzed to determine phase angle and amplitude. The theo- 
retical advection was determined by taking the mean zonal 
component of the wind (u) over a given time period. A per- 
centage of the actual advection/theoretical advection vs 
latitude for a 6-hour period is shown in figure io. 

Experiment IV. This experiment was designed to analyze 
the 8 effect as given in equation (4.2). It can be shown 
that the B effect is represented by the last term in equa- 
tion (4.2). The mean flow, B in equation (4.2), was set to 
zero which eliminates the large N-S pressure variation. The 
B effect is most pronounced for low wave numbers. Wave num- 
bers 4, 8 and 12 were run. According to equation (4.2) with 
B set to zero, wave number 4 should retrogress at 24 degrees 
per day, wave number 8 should retrogress at 8 degrees per . 
day and wave number 12 should retrogress at 4 degrees per 
day. The observed retrogressions are shown in figures 16, 
17 and 18. The actual retrogressions were somewhat less 


than expected due to the mean divergence which reduces the 
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magnitude of 8. The analyzed and 48-hour forecast surface 
pressure fields are shown in charts A through F. 

Experiment V. Wave numbers 4, 8 and 12 were run with 
the mean flow adjusted to give a phase speed of 10 degrees 
per day. As the wave number increases the second term in 
equation (4.2) has less contribution and the phase speed ap- 
proaches the initialized phase speed. Figures 19, 21 and 23 
show the observed phase angle vs latitude for wave numbers 
4, 8 and 12 respectively. As expected, the phase speed in- 
creased as wave number increased. The wave amplitude vs 
latitude are shown in figures 20, 22 and 24. In all three 
cases the wave form remained essentially constant over the 
48-hour forecast period. The analyzed and forecast surface 


eu 


pressurc fields are shown in charts G through L. 
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Phase angle (degrees longitude) vs latitude, for 
l6-point N-5 grid, wave number 4, phase speed 
u day and A = 7.0 x 10°. (Latitudes with zero 
wave amplitude are not included and time is 


given in hours.) 
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Figure 12. 


Phase angle (degrees longitude) vs latitude for 
16-point N-S grid, wave number 8, phase speed 
10° /day and A = 1.6 x 10°. (Latitudes with zero 
wave amplitude are not included and time is | 


given in hours.) 
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Figure 13. Phase angle (degrees longitude) vs latitude for 
16-point N-S grid, wave number 12, phase speed 
Os day and Alf 7.0 <x 10’. (Latitudes with zero 
wave amplitude are not included and time is 


given in hours. ) 
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Figure 15. Percentage of actual advection/theoretical advection vs 
latitude for wave number 4 with A = 7.0 x 10’, phase speed 
10°/day, Coriolis = 0, pressure gradient = 0, vertical 


advection = O and terrain pressure and temperature constant. 
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Figure 16. 
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Phase angle (degrees longitude) vs latitude for 46-point N-S 
grid, wave number 4, phase speed -24°/day, B=0 (no mean 
flow) and A = 7.0 x 10 (Latitudes with zero amplitude are 


not included and time is given in hours.) 
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Chart A. Initial surface pressure analysis for 46-point N-S grid, wave 
number 4, phase speed -24?/day, B=0 (no mean fiow) and 


A = 7.0 x 10?. 
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A Loss 
Chart B. 48-hour surface pressure forecast for 46-point N-S grid, wave 


number 4, phase speed -24° /day, B=0 (no mean flow) and 
A - 7.0 x 107. 
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Figure 17. Phase angle (degrees longitude) vs latitude for 46-point N-S 
grid, wave number 8, phase speed -8°/day, B=0 (no mean flow) 
and A = 1.6 * 10°. (Latitudes with zero amplitude are not 


included and time is given in hours.) 
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Chart C. Initial surface pressure analysis for 46-point N-S erid, 
wave number 8, phase speed -8°/day, B=0 (no mean flow) and 


A = 1.6 x 10°. 


60 





790'N 





_.-90'S 
Chart D. 48-hour surface pressure forecast for 46-point N-S grid, 


wave number 8, phase speed -8°/day, B=0 (no mean flow) and 


Aw 16 x 10°. 
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Figure 18. Phase angle (degrees longitude) vs latitude for 46-point N-S 
grid, wave number 12, phase speed -4°/day, B=0 (no mean flow) 
and A= 7.0 x 10’. (Latitudes with zero amplitude are not 


included and time is given in hours.) 
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wave number 12, phase speed -4^/day, B-0 (no mean flow) 


AO x 10^. 


63 





1002 


1006 






1002 






1002 
! 
| 
x 


Chart F. 48-hour surface pressure forecast for 46-point N-S grid, 


wave number 12, phase speed -4°/day, B=0 (no mean flow) and 


90°S 


A 0505506107. 
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Figure 19. Phase angle (degrees longitude) vs latitude for 46-point N-S 
grid, wave number 4, phase speed 10°/day and A = 7.0 x 10’. 
(Latitudes with zero amplitude are not included and time is 


given in hours.) 
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Figure 20. Terrain pressure amplitude vs latitude for initial field 
and 48-hour forecast, wave number 4, phase speed 10°/day 


and A = 7.0 x 10’. 
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Chart G. Initial surface pressure analysis for 46-point N-S grid, wave 


number 4, phase speed 10°/day and A = 7.0 x 10°. 
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Chart H. 48-hour surface pressure forecast for 46-point N-S grid, wave 


number 4, phase speed 10°/day and A = 7.0 x 10’. 
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Phase angle (degrees longitude) vs latitude for 46-point N-S 
grid, wave number 8, phase speed 10°/day and A = 1.6 x 10°. 
(Latitudes with zero amplitude are not included and time is 


given in hours.) 
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Figure 22. Terrain pressure amplitude vs latitude for initial field 


and 48-hour forecast, wave number 8, phase speed 10°/day 


and A = 1.6 x 10°. 
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Chart I. Initial surface pressure analysis for 46-point N-S grid, 


wave number 8, phase speed 10°/day and A = 1.6 x 10°. 
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Chart J. 48-hour surface pressure forecast for 46-point N-S grid, 


wave number 8, phase speed 10°/day and A = 1.6 x 10°. 
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Figure 23. Phase angle (degrees longitude) vs latitude for 46-point N-S 
grid, wave number 12, phase speed 10°/day and A = 7.0 x 10’. 


(Latitudes with zero amplitude are not included and time is 


given in hours.) 
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Figure 24. Terrain pressure amplitude vs latitude for initial field 
and 48-hour forecast, wave number 12, phase speed 10°/day 


and A = 7.0 x 10’. 
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Chart K. Initial surface pressure analysis for 46-point N-S grid, 


7 
wave number 12, phase speed 10°/day and A = 7.0 x 10’. 
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Chart L. 48-hour surface pressure forecast for 46-point N-S grid, 


wave number 12, phase speed 10°/day and A = 7.0 X 1o 
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VI. CONCLUSION 


The model remained well-behaved throughout all experi- 
ments and the phase speeds were less than for the non-diver- 
gent model, however, they equaled or exceeded theoretical 
expectations as described by Williams (1972). The model, as 
it is presently constructed, appears to have the desired 
horizontal flexibility as evidenced by the fact that the ex- 
periments were carried out with different grid spacings. 

The vertical flexibility, although built in, was not tested 
in this report. The simplified nature of the model pre- 
cluded direct comparison with existing global models such as 
FNWC's 5-level model. 

Future research efforts with this model should most 
likely include expansion to 5-levels, treatment of cross 


polar flow and initialization with real data. 
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